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Abstract 



We propose a new method for the problems of computing free en- 
ergy and surface pressure for various statistical mechanics models on a 
lattice Z''. Our method is based on representing the free energy and 
surface pressure in terms of certain marginal probabilities in a suitably 
tS^ modified sublattice of Z"^. Then recent deterministic algorithms for com- 

rrt puting marginal probabilities are used to obtain numerical estimates of 

a the quantities of interest. The method works under the assumption of 

I I Strong Spatial Mixing (SSP), which is a form of a correlation decay. 

We illustrate our method for the hard-core and monomer-dimer mod- 
T-H els, and improve several earlier estimates. For example we show that the 

^ exponent of the monomer-dimer coverings of Z^ belongs to the interval 

^~~^ [0.78595,0.78599], improving best previously known estimate of (approx- 

•^ imately) [0.7850,0.7862] obtained in [FP05] . [FKLM] . Moreover, we show 

^•^ that given a target additive error e > 0, the computational effort of our 

. method for these two models is (1/e)^^^ both for free energy and surface 

f pressure. In contrast, prior methods, such as transfer matrix method, 

^^ require exp ((1/e)'^'^') computation effort. 
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1 Introduction 

It is a classical fact in statistical physics that the logarithm of the partition func- 
tion of a general statistical mechanics model on \—'n^ri\'^ C Z'^, appropriately 
rescaled, has a well-defined limit as n ^- oo |Geo88| . |Sim93] . This limit is called 
free energy or pressure (the difference is in normalizing constant). Similarly, the 
limit of its first order correction, surface pressure, exists, though such limit de- 
pends on the shape of the underlying finite box and the existing proofs assume 
soft core interactions, H < oo (see Section [2| [ Sim93| . It is a different matter 
to compute these limits, and this question interested researchers both in the 
statistical physics and combinatorics communities. In very special cases, free 
energy can be computed analytically. The most widely known example is Fisher- 
Kasteleyn-Temperley's formula for dimer model on Z^ |Fis61] . |Kas61| . |TF61] . 
This formula was used in |She06j for obtaining complete parametrization of 
different Gibbs measures for the pure dimer model on Z^. The formula, un- 
fortunately, extends neither to a dimer model in other dimensions nor to the 
closely related monomer-dimer model. Another example of an exactly solvable 
model is hard-core model on a hexagonal lattice |Bax80] . 

Short of these special cases, the existing methods for computing free energy 
mostly rely on numerical approximations. These include randomized methods 
such as Monte-Carlo [JS97 , [KRS96 and deterministic methods such as transfer 
matrices. The Monte-Carlo method can be used to estimate free energy in finite 
graphs, say [—n, n]'^, with some probabilistic approximation guarantee, provided 
that some underlying Markov chain is rapidly mixing. One then has to relate 
finite graph to infinite lattices to approximate free energy for an infinite lattice. 
This can be done since the rate of convergence in the definition of the free energy 
is known to be 0(l/n) different from the log-partition function on [— n, n]'*, and 
the constant in O(-) can be explicitly bounded. The drawback of this method 
is its dependence on the sampling error. Transfer matrix method, on the other 
hand is a deterministic method and provides rigorous bounds on the free energy. 



It is based on considering an infinite strip [—Ti,n]'^~^ x Z and then identifying 
the number of different ways two configurations on [— n, n]''"^ can match. One 
then constructs an [—n,n]'^~^ by [—n,n]'^~^ matrix and the spectral radius of 
this transfer matrix can be related to the growth rate of the partition function 
on [—n,n]'^~^ x [—N,N] as a function of N, namely the partition function on 
[— n,n]''~^ X Z. Then by making n sufficiently large, bounds on the free energy 
on the entire lattice Z'^ can be obtained. The construction of such transfer 
matrix requires time exp(0(n''~^)) time. Since the convergence rate of the free 
energy to its limit wrt n is 0{l/n), then in order to get a target additive error e, 
the transfer matrix method requires time exp(0((l/e)'^~-'^). Namely, this is time 
required to construct an e-length interval containing the actual value of the free 
energy. While, for low dimensions d, this is a substantial saving over a brute 
force method, which would simply compute free energy on [—n,n]'^ for n > 1/e 
and requires time exp(0((l/e)''), the method stops being effective for larger 
d. Some additional computational savings can be achieved using underlying 
automorphisms group structure, see Friedland and Peled |FP05j . but we are 
not aware of any formal analysis of the computational savings produced by 
this method. The complexity of computing surface pressure using th e tr ansfer 



matrix method further increases to exp(0((l/e)^'^~^), see Subsection 3.3 

The transfer matrix method was used to obtain some of the best known 
bounds. For the hard-core model (see Section l4| with activity A = 1, which 
corresponds to counting independent sets in Z"*, the exponent of the free energy 
is known to be in the range [1.503047782,1.5035148], as obtained by Calkin 
and Wilf |CW98| . A far more accurate but non-rigorous estimate was ob- 
tained by Baxter |Bax99j . Similarly, the transfer matrix method was used to 
compute the free energy of the monomer-dimer (matchings) model (see Sec- 
tion [5]). The problem has a long history. Earlier studies include Hammersley 
|Ham66a] . |Ham66b] . Hammersley and Menon [HM70] . Baxter Bax68J, where 
some non-rigorous estimates and crude bounds were obtained. Recently Fried- 
land and Peled [FP05] obtained rigorously a range 0.6627989727± 0.0000000001 
for d = 2 and [0.7653,0.7863] for d= 3 using the transfer matrix method. The 
lower bound was later tightened to 0.7845 using the Friedland-Tveberg inequal- 
ity |FG08| . This inequality provides a bound for general regular graphs. This 
bound, while quite accurate for the case Z^, is not improvable by running some 
numerical procedure longer or on a faster machine. A tighter lower bound 
0.7849602275 was recently obtained using techniques related to the asymptotic 
matching conjectures FKLM . A bit earlier a similar non-rigorous estimate 
[0.7833, 0.7861] was obtained in ,HLLB08j . by reduction from a permanent prob- 
lem. Exact values for monomer-dimer entropy on two dimensional infinite strip 
were obtained by Kong [Kon07J . We are not aware of any computational esti- 
mates of surface pressure for these or of any related statistical mechanics models. 
In this paper we propose a completely new approach for the problem of 
computing numerically free energy and surface pressure. Our approach takes 
advantage of the fact that some of these models, including the two models above, 
are in the so-called uniqueness regime. Namely, the Gibbs measure on the infi- 
nite lattice is unique. This is an implication of the Strong Spatial Mixing (SSM) 
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Figure 1: G^ is obtained by removing all dark nodes 



property [BMP99], see Subsection 2.4 for the definition. This property asserts 
that the marginal probability that a node v attains a particular spin value ay 
is asymptotically independent from the spin configurations for nodes u which 
are far away from v. Our main theoretical result is Theorem IT] which provides a 
surprisingly simple representation of the free energy in terms of such a marginal 
probability, provided that the SSM holds. This representation is particularly 
easy to explain for the special case of hard-core model in two dimensions. Con- 
sider the subset Z^q of Z^ consisting of points (wi,i'2) with either ti2 < or 
V2 = 0, wi < 0. Our representation theorem states that the free energy equals 
log(l/p*) where p* is the probability that a random independent set in Z^q con- 
tains the origin (0, 0). As the model satisfies SSM (see SectionB this probability 
is well-defined. The idea of the proof is simple. Let vi ~i V2 -< ■ ■ ■ -< i'(2n+i)2 
be the lexicographic ordering of nodes in [— n,n]^. That is (ui,W2) >- (mi,M2) 
iff V2 > U2 or V2 — U2 and vi > ui. Observe the following telescoping identity 

Z^^ — Yu=o ~ ■^(*)/^(* + 1); where Z{i) is the number of independents 
sets in the subgraph Gi of [— n,n]^ when all nodes Vj,j > i and incident edges 
are removed, see Figure IT] The convention Z{0) = 1 is used. Observe that 
Z{i)/Z{i + 1) is the probability that a randomly chosen independent set in G^ 
does not contain Wi+i. This probability is then approximately p* for "most" 
of the nodes in [— n,n]^. The required representation is then obtained by tak- 
ing logarithms of both sides. Such representation is called cavity method in 
statistical physics and has been used heavily for analyzing statistical models on 
sparse random graphs, !Ald01'. [RBMM04j . |MP03] . [MP87] . where one can write 
certain recursive distributional equations satisfied by cavity values. Here due 
to a particular sequencing of removed vertices, we call our approach sequential 
cavity method. 

We establish a similar representation result for surface pressure for rectan- 



gular shapes again under the assumption of (exponential) SSM. We note that 
our results provide an independent proof of the existence of the free energy 
and surface pressure, not relying on the sub-additivity arguments, albeit in the 
special case of SSM. Theorem [T] reduces the problem of computing free energy 
and surface pressure to the one of computing marginal probabilities, and this is 
done using recent deterministic algorithms for computing such marginals in cer- 
tain models satisfying SSM property [Wei06]. [GK57] .[BGK+0r . These results 
are based on establishing even stronger property, namely correlation decay on a 
computation tree and lead to efficient algorithms for computing such marginals. 
The original intention of these papers was constructing scalable (polynomial 
time) approximation algorithms for computing partition functions on general fi- 
nite graphs, but an intermediate problem of computing approximately marginal 
probabilities is solved. We have implemented these algorithm for the hard- 
core |Wei06] and monomer-dimer BGK+0'^ models in the special case of U^. 



Using our approach we improve existing bounds for these models. For example 
we show that the exponent of the free-energy of the hard-core model in two di- 
mensions is in the range [1.503034, 1.503058]. While our lower bound is weaker 
than the previous best known bound 1.503047782 |CW98j (see above), which is 
already quite close to a believed estimate |Bax99| . our upper bound improves the 
earlier best bound 1.5035148 |CW98j . For the case of monomer-dimer model 
in three dimensions we obtained a range [0.78595,0.78599], substantially im- 
proving earlier bounds |FP05] . [FG08] . |FKLM| (see above). Further we show 
that the numerical complexity of obtaining e-additive approximation using our 
approach is (l/e)*^^^^ for both free energy and surface pressure. The constant 
in 0(1) may depend on model parameters and dimensions. For example, for 
the case of monomer-dimer model, this constant is Cd/i logd, where C is some 
universal constant (Proposition El . This is a substantial improvement over the 
computation efforts exp(0((l/e) ) and exp(0((l/e)^'^~^) of the transfer ma- 
trix method for computing free energy and surface pressure. 

The rest of the paper is organized as follows. In the following section we pro- 
vide a necessary background on Gibbs measures on general graphs and lattices, 
free energy, surface pressure, and define an important notion - (exponential) 
Strong Spatial Mixing. Our main theoretical result is Theorem IT] which repre- 
sents free energy and surface pressure in terms of marginal probabilities. This 
result and its several variations are stated and proven in Section [3] Sections [4] 
and |5] are devoted to application of Theorem [l] and its variations to the prob- 
lem of numerically estimating free energy for hard-core and monomer-dimer 
models specifically. Additionally, in these sections we compare the algorithmic 
complexity of our method with the complexity of the transfer matrix method. 
Concluding thoughts and open questions are in Section [6] 



2 Model, assumptions and notations 

2.1 Finite and locally finite graphs 

Consider a finite or infinite locally finite simple undirected graph G with node 
set V and edge set E C V xV. The locally-finite property means every node is 
connected to only finitely many neighbors. The graph is undirected and simple 
(no loops, no multiple edges). For every v,u Cz V let d{u,v) be the length 
(number of edges) in the shortest path connecting u and v. For every node v, 
we let N{v) stand for the set of neighbors of v. N{v) = {u : {u, v) £ E}. We will 
write Nq^v) when we need to emphasize the underlying graph G. The quantity 
A = Ag = maxt, |Af(w)| is called the degree of the graph. For every r E Z+ 
let Briv) = {u : d{v,u) < r}. For every set A d V, let Br{A) = U^(=ABr{v). 
Thus locally finite property means |i?r(v)| < oo for all v,r. Given A C V, let 
dA ^ {v e A : nIv) n A= 7^ 0} and let drA = d{Br{A)). For every subset 
A d V, wc have an induced subgraph obtained by taking nodes in A and all 
edges {v,u) G E A^. The corresponding edge set of the induced graph is 
denoted by E{A). 

Let K+ (Z+) denote the set of all non-negative real (integer) values. Le 
M>o(Z>o) be the set of all positive real (integer) values. Our main example 
of an infinite locally-finite graph is the d-dimensional lattice Z'^ with V = 
{{vi, . . . ,Vd) : Vi € Z and E = {{v,u) € V'^ : \\v — u\\ = 1}, where ||w|| = 
X^KiXd \^j\- We denote the origin (0, 0, ... , 0) by for short. Given a vector 
a = {ai, . . . ,ad) G M'J_, and v = {vi, . . . ,Vd) G Z'', let Ban{v) = {u £ 1/ : 
Wj — Vj\ < ttjU, j = 1,2, .. . ,d}. In the special case G = 'Z'^,v = 0, we write 
Br instead of Br{0) and Ban instead of i?a„(0). For each j < d,k £ Z+, let 
Z^j.^ = {v e Z'^ : Vj < k},Z'^,^_ ^ {v £ Z'^ : v^ > ~k}. Each of these 
is a d-dimensional half-plane in Z''. Let -< denote a lexicographic full order 
on Z''. Namely, u -< m iff either Ud > Vd or 3k £ {1,2, . . . ,d — 1} such that 
Uk > Vk and Vj = Uj,j = k + 1, . . . ,d. For every v £ Z"^ and fc = 1 , 2, . . . , d, let 
^i^v = {u £ Z"^ : u ~< v}L} {v}. Similarly we define Z^ j,^^^^ and Z^^, _^^^ with 
Z^^,^ _^ and Z'^^^^_ replacing Z'^. 

Throughout the paper we write /(n) == 0{g{n)) and f{n) ~ o{g{n)), n G Z+ 
if /(n) < Cg{n), respectively f{n)/g{n) -^ 0, for all n, for some constant C. 
This constant may in general depend on model parameters such as H, h (see the 
next section) or dimension d. However, in some places the constant is universal, 
namely independent from any model parameters. We will explicitly say so if 
this is the case. 

2.2 Gibbs measures 

Consider a finite set of spin values x = {^ij • • • i^q}, a Hamiltonian function 
H : x^ ^ ^^ {00}, and an external field ft, : % ^- M. Given a graph G = (V, E) 
we consider the associated spin configuration space f2 = x equipped with 
product cr- field !F. If the graph G is finite, a probability measure P on (fi, !F) is 



defined to be Gibbs measure if for every spin assignment (s„) S x^ 

PK = s„, Vwey) = Z-iexp(- Yl -ff(s.,s„)-^/i(s,)), (1) 

{v,u)eE vev 

where Z is the normahzing partition function: 

Z= ^ exp(- ^ H{s^,Su) -^hisy)). 

(siOex'^i {v,u)eE vev 

We will often write P^ and Zq in order to emphasize the underlying graph. The 
case H(a, a') = oo corresponds to a hard-core constraint prohibiting assigning 
a G X Etnd a' G X to neighbors. The possibility of such hard-core constraints 
is important for us when we discuss the problems of counting independent sets 
and matchings. 

When G is infinite, a probability measure P(-) is defined to be Gibbs measure, 
if it satisfies the following spatial Markovian property. For every finite A d V 
and every spin assignment (s„) € x on the boundary dA of A, the conditional 
probability measure P(-|(s^)) on the finite graph {V{A),E{A)) induced by A is 
a Gibbs measure with the same H as the original graph and external field h' 
given as K{a) = h{a) + Z]«eA'(i>)\A '*«('5") foi' every a £ x- Here every node 
u € N{v) \ A belongs to dA and thus its spin value s„ is well defined. This is 
called spatial Markovian property of Gibbs measure. One can construct Gibbs 
measures as a weak limit of Gibbs measures on cylinder sets in !F obtained 
from finite induced subgraphs of G, see |Sim93] . rGeo88] for details. Generally 
there are multiple Gibbs measures and the space of Gibbs measures of G is 
denoted by M or Mg- We say that the model {G,H,h) is in the uniqueness 
regime if A4 consists of a unique measure A4 = {¥}. Most of the results in 
this paper correspond to the uniqueness case, and more specifically to the case 
of Strong Spatial Mixing (SSM) defined below. Given a model (G, H, h) and a 
subset ^ C T^, we have a naturally defined induced submodel on the the induced 
subgraph G{A) = {V{A),E{A)), given by the same x, -ff and h (but different 
partition function). In order to emphasize the underlying subgraph G we write, 
with abuse of notation, Za and P^ for the partition function and the Gibbs 
measure on the subsystem {G{A),H,h), when it is unique. By default we drop 
the subscript when the underlying graph is the entire lattice Z'^. 

The definitions above is "node" centered: the spins are associated with nodes 
of a graph. In order to study the monomer-dimer model we need to consider a 
similar "edge" model where spins are associated with edges. Thus given a finite 
graph G = {V,E), and a finite set of spin values x. we consider a probability 
space fl — x''^'- Then ([l| is restated as follows. Write e ^ e' if edges e, e' are 
distinct and incident (share a node). For every spin assignment (s^) G x we 
assign probability measure 

P((Te-Se, Vee^) = Z-iexp(- J2 H{s,,s,,)-Y,Hse)), (2) 



It is clear that "edge" model can be reduced to "node" model by considering 
a line graph of G: the nodes of this graph are edges of G and two nodes e, e' 
form an edge in the line graph if and only if e ~ e'. For simplicity, however, we 
prefer not to switch to the line graph model. 

2.3 Free energy, pressure and surface pressure 

Given G — Z'^, Hamiltonian H and an external field h, consider an arbitrary 
infinite sequence of finite subsets S Ai C A2 C • • • C Z'', such that the 
sequence r„ = max{r : B^ C A^} diverges to infinity as n — > 00. Consider the 
corresponding sequence of Gibbs measures Pa„ and partition functions Za„ on 
the graphs induced by A„. It follows from sub-additivity property of partition 
functions that the limit 

P(d,i/,/.)^ lim i^^ (3) 

n^oo |A„| 

exists and is independent from the choice of the sequence of subsets [Sim93', 
|Geo88j . This quantity is called pressure. Given a positive /3 > one usually 
also considers limits lim„^oo log.^A„/(/3|A„|) where H and h are replaced by 
PH and /3/i, respectively. The corresponding limit is called free energy. For our 
purposes, this difference is insubstantial, as it is just a matter of redefining H, h 
and changing the normalization. We will mostly use the term free energy. 

The first order correction to the limit dsl is known as surface pressure and is 
defined as follows. Unlike pressure, this quantity is "shape" dependent. In this 
paper we only consider the case of the rectangular shape, although our results 
can extended to more complicated shapes as well. Given a vector a = (ai) € 
Rio, let A,{a) ^ nfc^,(2afc), J = 1,2,..., d, and let A{a) = Ej<dMa)- 
Observe that surface area of Ban (the number of boundary nodes in Ban) is 
A{a)n'^~^ + o{n'^~^). The surface pressure is defined as 

sV{d,H,h,a) = 

= Ym,^A-\a)n-''+ A log ZB^^io)-n''V{d,H,h) \{ {2ai)\ . (4) 

"~*°° \ l<j<d J 

The surface pressure is interpreted as the first order correction of the free energy, 
scaled by the surface area. The existence of this limit is a classical fact |Sim93] 
for the case maxg^j/g^ H{s, s') < 00. The case hard-core case i? = 00 is trickier 
and we are not aware of results on the existence of the limit for this case. Our 
proofs however imply the existence of the limit under the assumption of SSM, 
see below. 

2.4 Strong Spatial Mixing 

The following technical assumption is needed for our analysis. The assumption 
essentially says that there is always a choice of a spin value s* E x which allows 



any other choice of spin values for its neighbors without making the Hamiltonian 
infinite. 

Assumption 1. There exist s* d x such that laaxsex I^(*j**)I < °°- 
The assumption imphes that under any Gibbs measure, for every node v 

c* = minP(cr„ = s*|cr„ = s„,Vu G N{v)) > (5) 

where the minimum is over all possible spin assignments (s„) € t^I"(")I of the 
neighbors of v. We now introduce the Strong Spatial Mixing assumption, which, 
in particular, implies the uniqueness of the Gibbs measure on an infinite graph. 
We first informally discuss this notion. Loosely speaking, a model (G, H, h) 
exhibits just spatial mixing if for every finite set A d G the joint probability 
law for spins in A is asymptotically independent from the values of spins which 
are far away from A. A form of this condition is known to be equivalent to the 
uniqueness of the Gibbs measure. Instead, a strong spatial mixing means the 
same except that values for any other subset _B C Z'' are allowed to be fixed. In 
other words, strong spatial mixing is spatial mixing, but for the reduced model 
on y \ _B when values of spins in B are fixed and values of the external field h on 
the remaining graph are appropriately modified. Strong spatial mixing is strictly 
stronger than spatial mixing (hence the separate definition), see |BMP99| for 
counterexamples and discussion. 

We now introduce the definition formally. 

Definition 1. (G, H, h) satisfies the Strong Spatial Mixing property if for every 
finite set X C V^(G), possibly infinite set Y e V{G) there exists a function R{r) 
satisfying lim,.^oo Ri^) = 0, such that for every positive integer r 

{(Ty = sl,v e X\ay ^ sl,v eY;cry = sl,v e drX) 

-P(ct, = sl.v e X\ay = sl,ve Y;<Jy = st,v e a,X)| < R{r), 

where the maximum is over all possible spin assignments (s^) € X i (s^) G 
x'^', (s^), (s^) € x' ■ (G, -ff, /i) satisfies exponential strong spatial mixing if 
there exist k,7 > such that R{r) < Kexp(—jr) for all r >0. 

It is known that if (G, H, h) exhibits SSM, then the Gibbs measure is unique: 
M = {V}, [BMP99^ . 

3 Sequential cavity method 

In this section we present our main theoretical result: the representation of the 
free energy and surface pressure on Z'' in terms of some conditional marginal 
probabilities P(cro = s\-) defined on suitably modified subsets of Z*^, 'Zj^k,+ ,'^j,k,-- 
The idea is to sequentially remove nodes from a rectangle Ban one by one and ob- 
serve that the log-partition function can be written in terms of the log-marginal 
probabilities (cavity) of the removed nodes. 
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3.1 Representation theorem 

Given a graph G = (VtE) with some full order >-, a node v and s € Xi lot £„^s 
denote the event "a^ — s for all u '^ v,u € V" . By convention we assume that 
£v^s is the full event fi, if the set u )^ v,u G F is empty. In the special case when 
G = Z"*, w = and >~ is the lexicographic order, this event is denoted by Eg, see 
Figure [2] for the case d = 2. We remind the reader, that we drop the subscripts in 
Pg, Zq when the underlying graph is Z'^ and SSM holds guaranteeing uniqueness 
of the Gibbs measure. 

Theorem 1. Suppose {!/'■. H.h) satisfies the SSM property and the Assump- 
tionUl holds. Then 

V{d,H,h) = -logP((To = s*\£s^)-dH{s*,s*) - h{s*). (6) 

//, in addition, SSM is exponential, then for every a 

sV{d,H,h,a) = -H{s\s*)+ 

. / -. OQ 

-logP^. (ao = s*|£,.) 

Moreover, the infinite sum is convergent. 

Remark : We stress that our proof docs not rely on the existence of a limit 
in (|3| and thus provides an independent proof for it, albeit in the restricted case 
of SSM. Preempting the formal discussion, the convergence of the infinite sum 
in drl is an almost immediate consequence of exponential SSM, which gives 

|Pz. (ao = s*|£s-)-PK = s*|£,.)l <Kexp(-7fc). 

A similar conclusion for difference of the logarithms is obtained then using 
Assumption [Tl This observation implies that e-additive approximation of this 
infinite sum is obtained by considering the partial sum with first 0(log(l/e)) 
terms. As we shall discuss in the following subsection, the main import from this 
observation is that the overhead of computing surface pressure vs free energy is 
very small. 

Proof. We first prove (IbI. Let vi)~ V2> . ■ . >- ^|s„| be the labeling of nodes in 
Bn according to the lexicographic order. Note that the number of edges in i?„ 
is d|S„|+o(|B„|). We have 

Pi3„(fT, = s*,VweB„) 

= Zgle^^{-d\Bn\H{s*,s*)^o{\B^\)H{s*,s*)^\B,,\h{s*)), 

from which we infer that 

Zgl = e^^{d\Br,\H{s\s*)+o{\Br,\)H{s\s*) + |B„|/i(s*))x 

xPs„(a„ = s*,VveB„). (8) 
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Figure 2: Event Eg* on Z^. Every dark node is assigned spin s* 

On the other hand, by telescoping property 

Ps„K = s*,Vf eB„)= '[[PBA'^v = s*\au^s*,yuyv,uGBn) (9) 



veB„ 

n 

t;6B„ 



,(cr„ = S*\£y^s*) 



Fix e > and find r = r(e) such that according to Definition fl] R{r) < e 
for X = {0},Y = {u y 0}. Let B° ^ {v G B„ : Br{v) C B„}. Observe that 
|i3°|/|i?„| — ^ 1 as n ^ oo. By the choice of r we have for every v G B° 



(ay = s* \£y^s' ) - Pb„ (cr„ = s* |£^,s- ) < e. 
By translation invariance we have 

¥{ay = s*\£y^s') ^n<7o = s*\£s') 
For every w G _B,i \ 5° we have the generic lower bound (Is]) 



(10) 



(11) 



which is strictly positive by Assumption [Tl A similar inequality with the same 
constant holds for P(cro — s*\£s')- Combining this with (Is]) and (|9| we obtain 

Zsl = exp{d\Bn\H{s*,s*)+o{\Bn\)H{s*,s*) + |S„|/i(s*))x 

veB° veBnXB" 

> exp {d\Bn\H{s*,s*) + o{\B^\)H{s*,s*) + |S„|/i(s*)) 

(P(ao = .1£.-)-^)"'°'(c*)l^"\^°l 
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Since |i?°|/|i3„| ^ 1 as n — + cx3 and c* > 0, then we obtain 



lini inf ■ 



IB. 



'f- > dH{3*,s*) + h{s*) + log(P(ao - s*\£s') - e) 



Recalling P(cro = s*\£s<.) > c* > 0, since e was arbitrary, we conclude 
lim inf ■ 



■ - °^ ^" >H{s*,s*)d+h{s*)+logP{ao = s*\£s'). 



B„ 



Similarly we show 



lim sup 



logZ^ 

l^nl 



< H{s*,s*)d+ h{s*)+ log¥{ao = s*\£,,), 



where for the case v E Bn \ B°^ we use a trivial inequality Ps„ {(Jy — s* \£v,s* ) < 1 
in place of (111. We obtain 



lim ^^ = H{s*,s*)d+h{s*)+\ogV{aa = s*\£,.), (12) 

This concludes the proof of (|6| . 

Now we estabhsh Q. Thus consider a rectangle Ban- The proof is based 
on a more refined estimates for the elements of the telescoping product (|9l with 
Ban replacing i3„ . We begin by refinement of ([S]) . The number of edges in Ban 
is (i|i3a„| — A{a)n'^~^ + o{n''-~^). Repeating the derivation of (|8| and (pi), we 
obtain 



exp {d\Ban\H{s\s*) - A{a)n''-^H{s*,s*) + h{s*)\Ban\ + o(n'^"i)) 



Yl Pb_((7, = s* !£,,,.) 



(13) 



veBa 



Let 



B°„ = jt; e Ban ■ \ain~v,\,\ain + v,\ > Clogn, 1 <i <dj, 
and for each j — I, . . . ,d and < k < C log n let 
Ban.j,k,+ = Iv e Ban ■ v-j = [flj-TiJ - fc; \ain - v^\,\ain + v.^ > Clogn, i 7^ j|, 



B, 



a7i,j,k.- 



{« 



^iveBa 



'lajn\ + k; \ain — Vi\, jajTi + Vi\ > Clogn, i ^ j>- 



Here C > is a large yet unspecified constant. In other words, Ban,j.k,+ U 
Ban,j.k.- consists of nodes which have distance k from from the boundary of 
Ban in the j coordinate, but have distance at least Clogn in all the other 
coordinates i ^ j. Observe that the sets Ban,j,k.+ , Ban.j,k,- are non-intersecting, 
for sufficiently large n. Observe also that 



Ban \ [Ban U ^l<j<d,0<k<C \ogniBan,j,k,+ U Ban,j,k~)) 

= 0{n'^-^\og^n). 



(14) 
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-ain 



Bn 



V2 
02" 



> Clogn 



(0,0) 



-02" 



->o<-»-— •---»--♦-> • 

/^ > Clogn 



ain 



->Dl 



> Clogn 



Figure 3: Event Eg* on'Z'^ , ^i , . Every dark node is assigned spin s* 



For every v E _B°„ we have by exponential SSM 

|P(a„ = s*\£^,s')-Vb,A<^v - s*\£,^s')\ < «:exp(-7Clogn) = 0(n-''^). 
By translation invariance we have 

PK -s* !£„,,. )-iPK-s*|fs-) 

Recalling P(iTo = s*\£s') > c*, we obtain using the Taylor expansion, 

I logPB„.(fT„ = s*|£,,,.) - logP(ao = s*\£,,)\ = O(n-T^). (15) 

Fix V G Ban.j,k,+ - Since Vj = [fljnj — k and for every i =/^ j , Vi is at least 
Clogn away from the boundary of Ban, then by the exponential SSM property 
we have 

^{ay = s*\£y,s-*) -Vj,d {(7y = s*\£.u^s') < K cxp ( " 7C log Ti) = O (n"^"^ ) , 

(see Figure [3] for the two-dimensional illustration with j — 2). Recalling Vj = 
[fljnj — fc and using translation invariance we have 

P^d {<J„ =S*\£,,s-*)= Pz^ {<JQ = S*\£s* ) 

Using again the Taylor expansion, we obtain 

logPB„„(a, = 5*|£,.,0-logPz<* (ao = s*|£.0 =0(n-T^). (16) 
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Similarly, if v G Ban,j,k.-, then 



-B_(a„-s*|£,,,0-IPz^ {(Jv = s* ]£,,,>) <exp(-7Clogn) = 0(n-^c^). 

i.-La^iiJ,- 



By translation invariance and since Vj ~ — [a^nj + fc, then 



j,-L"j"J.- 

and again applying Taylor expansion 



(O-Q = S*\£s'), 



logPB_(a„ = s*|£,.3.)-logPz^ (ao-s*|f,.) = O(n'T^). 

j,-fi,- 



(17) 



We now take log of both sides of ( 13 1 and divide by n'^ ^ to obtain 



A{a)H{s*,s*) + '-^h{s*) + o(l) 

1 ^ logPB_K = S*|£„,,.) 



,-d+ 



De-B° 



"'+' E E E ^og¥s^Aa.,^s*\£^^,.) 

1< j<dk<C log nveBan.j^k, - 



n 



Here we use ( 14 1 and the fact Pb„„ (ctu — s* \£v,s' ) > c* > for all v. Applying 
PI 



^ log Pb_ (a, = s* !£,,,.) = 



v£B° 



= |i?°JlogP(ao = s*\£,.) + |i?:„|0(n-^^) 

= |E,„|logP(ao = s*|£,.) "logPK = s*|fs-)x 

E E (|Ba„,,A+| + |i?a„,jU-.+ |)+0(«'"'), 

l<j"<di:<Clogn 



where we used (14 1 in the last inequahty and assume that C > I/7. Using the 
established identity ^ we conclude 



\0gZB^„ + \Bar^\V{d,H,h) 



TV 



d-1 



-A{a)H{s*,s*)+ 



-d+l 



■""'"' E E E (logPi3„„(fT„=.S*|£„,,.)-logP(^0 = S*|£,.)) 

""'^' E E E ( log Pb„„ K = s* !£:„,,.)- log P(ao = s*|£,.)) 



l<j<dfc<Clognt)eSa„,j,fc_, 



0(1) 
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Further applying (16) and (17 1 and using C > I/7, we obtain the following 
expression 



~A{a)H{s*,s*)+ 

l<j<dk<Clogn 

l<j<dk<C\ogn 
+ 0(1) 

We have 

\Ban^j,k,+ \ = Aj{a)n''-^ + o(n'^-i), \Ban,3,k,-\ = Ajia)n''-' + 0(7/-!). 
The resulting expression is then 

-A{a)H{s*,s*)+ 

+ J2 ^j(«) Y. {^og¥^.^^Jao = s*\£,,,.)~log¥{ao^s*\£s')) 

l<j<d fc<Clogn 

+ J2 ^j(") E {logV^^_^_{ao = s*\£,,s*)-logV{ao^s*\£,.))+o{l) 

l<j<d k<Clogn 

Applying the exponential SSM property and the Taylor expansion, we have for 
every k> C log n 



logP^. (ao = s*|£,.,.) -logP(ao = s*\£s'))\ < Kexp{--fk) 
logP^. JfTo = s*|£,,,0 - logP(ao = s*\£,,))\ < Kexp(-7fc). 



This means that we can replace the sums X]fc<cioKn "^ith infinite sums X]fc>0' 
with a resulting error 0(exp(— 7Clogn) = 0{n~'^'~^). Dividing by A{a) we 
obtain the result. D 

3.2 Extensions and variations 

In this subsection we establish several variations of the first part of Theorem [Tl 
namely the representation ([6| of the free energy in terms of the marginal prob- 
ability. First we extend identity ([6| for the case when we do not necessarily 
have SSM, but instead have an upper or lower bound on marginal probability 
P(cro = s* \£s*)- In this case we obtain an analogue of ^ in the form of inequal- 
ities. These inequalities will be useful for obtaining numerical bounds on free 
energy for hard-core model in dimensions d = 3,4 in Section [4j 

Given r > consider an arbitrary spin assignment (su) G ^i"^^i which 
is consistent with event £5.. Namely, Su = s* for all w >- 0, u G dBr- Let 
Pmax{r) {pmin{f)) be the maximum (minimum) of PB^(cro = s*|£s.,(a„)) when 
we vary over all such assignments. 
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Figure 4: Event £s*,cvcn on Z^. Every dark node is assigned spin s* 



Corollary 1. For every (Z"*, H, h) and r > 

-^ogpynax{r) - dH{s*,s*) ~ h{s*) 
<V{d,H,h) 



< - logprmnir) - dH{s*,S*) - h{s*). 



(18) 



While the result holds for arbitrary r the quality of the bounds presumably 
improves with increasing r. We will see in Section H that in some cases Pmax(^) 
and Prnin('') are fairly close for large r even though the model is outside of 
provably exponential SSM regime. 

Proof. The proof is a minor variation of the proof of (l6| . Instead of estimate 
(fTol) we use Pmin{r) < PB„(cru = s*\£v,s') < Pmax(?') for every n> r. D 



Our second variation is a " chess-pattern" version of Theorem IT] where in 
representation ([9| we sequentially remove only vertices with even sum of coor- 
dinates. As it turns out this version provides substantial gains in computing 
numerical estimates of free energy both for hard-core and monomer-dimer mod- 
els, though we do not have theoretical explanation for this gain. 

Let I'eyen = {^' = (^ij • ■ • i ^d) <= Z'^ : X]i ^i ^^ cvcn}. Similarly define Z|^jj. 
Given a subgraph G = {V, E) of Z'^, let f„,s,even denote the event au = s,yu >- 
v,u gV n Zg^gjj. The special case u = is denoted by £s,cvcn- See Figure for 
this chess-pattern version of the event £3. . 

Theorem 2. Suppose {Z'^,H,h) satisfies the SSM property and the Assump- 
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tion\7\ holds. The 



V{d,H,h) = --logP(ao = S*\Es',even) 



^log(^exp(-2di?(s,s*)-/i(s))). (19) 



2 

Proof. Let i?„,cvcn = Bn n 1\,^ci\- Observe that every node v E Bn \ Bn.cvcn has 
only neighbors in -Bn^cvon- Let A(u) denote the degree of v in _B„. As a result, 

PB„K = s*,Vt;eB„,evc„) = ^s' n ^cxp(-A(t>)i?(s,s*)-/i(s)), 

t>es„\B„,evo„ 'sex 

from which we obtain 

^b'= n (5Iexp(-A(t;)i?(s,s*)-/i(s)))"Vs„K = s*,VweB„,evcn) 

■ueB„\s„,„vc„ 'sex 

On the other hand, by telescoping property 

Pb„(o-u = S*yv e Bn,cvcn) = Yl ^B^io^v = S*k„ = S*, V U ^ U, U € Bn,cvon) 

fS-Bn oven 

(20) 

I'e-Bn.ovon 

The remainder of the proof is very similar to the proof of Theorem [T] and details 
are omitted. Notice that for |_B„| — o(|i?„|) nodes in Bn, the degree A{v) = 2d. 



Then in place of ( 12 ) we obtain 

V logos' V \Bn\B, 



lim-y^ = - lim '^" ^7-"' log (^exp(-2rfg(.,.*) ~ h{s))) 

n \B„\ n^co \B.„\ V -^ — ' / 



sex 

lim %7^1ogP(ao = .1f...evc„) 
n^oo iy„ 



^log (^exp(-2dif(s,S*) - h{s))) + JlogP(ao = S*|£,-,eve„) 



2 

s6X 



n 



Now let us present a version of Theorem [T] for the model ([2]) where spins 
are assigned to edges rather than nodes. We need this for application to the 
monomer-dimer model. Let i?i,,s. edges denote the event cr„„, — s,Vu "^ v,w G 
N{u), and let ft,, s. edges, even dcnotc the event au,w — SjVu )^ v,u G Z^^^^-^,w G 
N{u). Let i?s, edges, i^is, edges, even dcuotc the Same events when v — 0. 



17 



Theorem 3. Consider a model (Z'*, H, h) given by (fM) with spins assigned to 
edges. Suppose {Z'^,H,h) satisfies the SSM property and the Assumptions] 
holds. Then 

V{d,H,h) = -logP(a(o,.) = S*yv e N{0)\£s^, edges) 

-d{2d-l)H{s*,s*)-dh{s*) (21) 

= --logP(cr(o,„) = S*,\fv e N{0)\£s*,edges,even) 

~d{2d~l)H{s*,s*)-dh{s*). (22) 

Remark : Contrast this result with Theorem [U There the extra term 
log(^gg •) appears since every odd node surrounded by even nodes with preset 
spin value s* still has x choices for the spin selection. For the edge version this 
is not the case: all edges are preselected to take spin values s* . 



Proof. We first prove (21 1. Let E{Bn) denote the edge set of i?„. Note that 
E{Bn) = d\Bn\ + o{\Bn\) and the the number of edges with 2d — 2 incident 
edges is also d\Bn\ + o(|i?„|). The number of pairs of incident edges is then 
d{2d-l)\Bn\+o{\Bn\). We have 

¥BA^e = s*,yeeE{B„)) 

= Z^lcxp ( - d{2d-l)\B,,\H{s*,s*) - o{\B,,\)H{s*,s*) - d\Bn\h{s*)). 

On the other hand, by telescoping property 

¥BA<^e^s*,yeeE{Bn)) 

= n ^Br,{'^(v,u) = s*,Vu e N{v)\a(^u,w) = s*, y uy v,w e N{u), u,w e Bn) 
veB„ 

= Y\. ^B„{cr{v,u) = 3*, yue N{v)\£y^s*, edges)- 
v(EB„ 



The remainder of the proof of (21 1 is similar to the one of Theorem [T] and is 

omitted. 

Turning to ( 22 1 , observe that every edge in Z"^ has exactly one end point in 



ZJ?y(,,j. Then we have again by the telescoping property 

¥BA<^e = s*,yeeE{Bn)) 

= Y\. ^B„{(^iv,u) = S*,yu € N{v)\(T^u,w) = S*, y Uy V,W G N{u), U e Bn,cvcn,W G -B„) 

= W Vb,M(v,u) = S*,yu& N{v)\£y^s>^cdgcs,c-wcn)- 

l'eS„.cvon 

The remainder of the proof of (21 1 is similar to the one of Theorem IT] The fact 
\Bn,even\l\Bn\ ~^ 1/2 && Ti ^ CO Icads to a factor 1/2 in (22l. D 
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3.3 Applications and numerical complexity 

Theorem [l] reduces the problem of computing V and sV to the problem of 
computing conditional marginal probabilities P(cro = s*\Eg*). This is certainly 
not the only way to represent free energy and surface pressure in terms of 
marginal probabilities. For example, consider a modified system (G, PH, (3h) on 
a finite graph G, and observe that 



vev (v,u)eE 

both expectations are taken wrt the Gibbs measure P. Thus knowing the 
marginal probabilities V{ay) and joint probabilities P(crt,, CTm), (w, u) £ E lets us 
recover log Zq approximately in principle. Unfortunately, this means we have 
to integrate the answers over (3 € [0,1], which in practice has to be approxi- 
mated by summation. In order then to guarantee the target level of accuracy, 
one would have to control the derivatives of marginal probabilities wrt f3. Addi- 
tionally, one would have to compute marginal probabilities for a whole range of 
/3, over which the integration takes place. The advantage of the representation 
(l6| and ^ is that it allows us computing free energy and surface pressure by 
computing only one marginal probability P(cro = s*\£s*)- In this paper we will 
compute these marginal probabilities approximately using recent deterministic 
algorithms for computing such marginal probabilities in general graphs, where 
appropriately defined computation tree (see following chapters) satisfies expo- 
nential SSM. As we will see below, when we have such property, our method 
provides an additive e approximation of free energy in time (1/e)*^'^', where the 
constant 0(1) may depend on the model parameters and d. Let us now compare 
this performance with the performance of the transfer matrix method. While 
we are not aware of any systematic numerical complexity analysis of the trans- 
fer matrix method, it can be deduced from the following considerations. The 
transfer matrix method is based on first computing the partition function on a 
strip [~n,n]'^^^ x Z. The latter is done by constructing certain Ixl^^"'*'"'^-' by 
|^|(2n+i) transfer matrix. For the hard-core case with A = 1 (see Section Wk 
the matrix is — 1 with 1 corresponding to allowed pair of neighboring configura- 
tions and corresponding to pairs of configurations which are not allowed. The 
spectral radius of the transfer matrix is then used to deduce the growth rate of 
the partition function, namely the free energy. Constructing such matrix takes 
time exp(0(n''~^)). Since the partition function on [— n,n]''~^ x Z converges to 
the one of Z*^ at the rate 0{l/n) [Sim93j . then in order to achieve an additive 
error e > 0, one needs exp(0((l/e)''~^)) computation efi^ort. 

Suppose one then wishes to use this method to approximate the surface 
pressure with an additive error e. What is the required numerical effort? We 
are not aware of applications of the transfer matrix method for computing the 
surface pressure. Thus the reasonable alternative approach is to 1) select a rect- 
angle an for some large value n, 2) Compute the partition function Zan in this 
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rectangle using perhaps the brute force method and 3) Compute the approxi- 
mation P of the free energy V using perhaps the transfer matrix method. The 
approximate surface pressure is then obtained from the fact that the conver- 
gence rate in the limit Q is 0{l/n) jSim93j . Let us show that this approach 
requires exp (0((l/e)^'' )) numerical effort in order to obtain e additive ap- 
proximation, regardless of how quickly one is able to compute Zan- Indeed, 
observe that this approach leads to an error n'^\P — V\/n'^^^ = n\P — V\. Thus 
for the target e additive error, we need to set n > 1/e. But this further re- 
quires that (l/e)|P - P| = 0(e), namely \P -V\^ 0{e^). Thus we need to 
achieve e^ additive error accuracy in estimating the free energy. This requires 
exp (0((l/e)^'^~^)) per our earher calculations, and the assertion is established. 
This is a stark contrast with complexity (l/e)"^*^^^ of the method proposed in 
this paper both for free energy and surface pressure for the cases of hard-core 
and monomer-dimer models. However, one should note that our method takes 
explicit advantage of the exponential SSM, while the transform matrix method 
does not rely on this assumption. 

4 Hard-core (independent set) model 

The hard-core lattice gas model, commonly known as independent set model in 
combinatorics, is given by x = {0, 1}, iJ(0, 0) = i?(0, 1) = H{1, 0) = 0, iJ(l, 1) = 
oo, /i(0) = 0, h{\) = f3 for some parameter (3. Case of interest is /3 < as it 
corresponds to Gibbs measure putting larger weight on larger cardinality inde- 
pendent set. Choosing s* = we obtain that Assumption [T] holds. 

It is common to set A — exp{—(3) > and let A be the parameter of the 
hard-core model. This parameter is usually called activity. Note that in terms 
of A, for every finite graph G 



^,l(. 



^V:<7^ = 1}\ 



(23) 



where the sum is over all spin configurations {<Ty) e {0, 1}^ such that aya^ = 
for all {v, u) e E. Equivalently, a subset of nodes / C V{G) is called an 
independent (also sometimes called a stable set) if for no edge (m, v) we have 



both u and v belong to /. Then we may rewrite ( 23 1 as 



= E^'" (24) 



where the summation is over all independent sets of G. The summation in ( 24 1 is 
sometimes called an independent set polynomial in the combinatorics literature. 
From now on we let I denote the random independent set selected according to 
a Gibbs measure (multiply defined if there are many Gibbs measures). In the 
case of finite graph G, for every independent set / we have 

Pg(/ = /) = Zc^aI^i. 
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The special case A = 1 corresponds to a uniform distribution on the set of all 
independent sets in G. We denote the free energy and surface pressure on Z'^ by 
V{d,X) and sV{d,\,a), respectively. In the special case A = 1, the free energy 
is also the entropy of the Gibbs distribution, since it is uniform. Conditioning 
on spins taking value s* = simplifies significantly in the context of hard-core 
model as the following proposition shows, a simple proof of which we include 
for completeness. 

Proposition 1. // a hard-core model on a graph G = (V, E) satisfies SSM 
for some A, then so does any subgraph of G. The same assertion applies to 
exponential SSM. Moreover, for every Wi,W2 C V{G) the following identity 
holds with respect to the unique Gibbs measure. 

Fg{v e i\Wi n / = 0, VK2 c /) = Pg(w e I), (25) 

for every v (^ V \ {Wi U -61(^^2)), where G is the subgraph induced by nodes in 
V\{WiUBi{W2)). 

Proof. Fix any positive integer r and consider any spin assignment (s„),w G 
dBr{v) \ {Wi U Bi{W2)). Extend this to a spin assignment to entire dBr{v) by 
setting Su — I (that is u e J), for m G dBr{v) n W2 and s„ = (that is u ^ I) 

for u e dBr{v) n iWi U {Bi{W2) \ 1^2)) • Call this spin assignment S. Applying 
spatial Markovian property we have 

PG(a. = 0|ct„ = 0,w e W^i,(t„ = 1,w e W2,S) 

^ PG(g. ^ 0, a„ ^ 0, M G Wi,(Ju = 1, M € VK2, \S) 
Pg(ct„ - 0, u e T4^i, a„ = 1, w e W2\S) 

where Ii is the set of independent sets in Br{v) such that v ^ 1,1 C] Wi = 
0, W2 C / and u € I n dBr{v) iff s„ = 1 (according to S). The set I2 is defined 
similarly, except the condition v ^ I is dropped. Then the ratio is equal to 

where J4 is the set of independent subsets of -Br-(O) \ (Wi U Bi(W2)) such that 
u € I D dBr{v) iff s„ = 1, and X3 is defined similarly, except in addition v ^ I 
for every 1^1^. We recognize this ratio as P^icTy = 0\S). We conclude 

PG(a„ = 0|a„ = 0,u e VKi,a„ ^l,ue W2,S) - Pj|(a, == 0\S). 

Now consider any induced subgraph Gi of G and any two sets Wi,W2 in Gi. 
Applying a similar argument we obtain 

Pgi (a, = 0\au = 0,ueWi,<Ju = 1, u G W2, S) 

= Pg(ct. = 0|(T„ = 0,ueWiU{V\ yi),a„ = l,u e 14^2, 5) 
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where Gi is induced by nodes Vi \ {WiUBi{W2))- By SSM property, the second 
quantity has a hmit as r — *■ cx) which is independent from assignment S on the 
boundary dBr{v). Therefore the same apphes to the first and third quantities. 
The first conclusion impHes that Gi satisfies SSM. The second conclusion gives 



(25) when applied to Gi = G. D 



In light of Propoposition [T] we obtain the following simplification of Theo- 
rems [T] and [2] in the hard-core case. Let 

Z^O.evcn = Ko U {" G Kid : -i u} , 

see Figure |4] 

Corollary 2. Suppose the hard-core model on Z'^ satisfies SSM for a given A. 
Then 

r{d,X) = -logF^.jOiI) (26) 

= -^ logP^.^ _ (0 il) + l log(l + A). (27) 

Thus we now focus on developing an algorithm for numerically estimating 



marginal probabilities appearing in ( 26 ) and ( 27 1 



4.1 Recursion, sequential cavity algorithm and correlation 
decay 

Let us now introduce a recursion satisfied by the hard-core model. This identity 
in a different form using a self-avoiding tree construction was established recently 
by Weitz |Wei06j . We repeat here some of the developments in |Wei06j , with 
some minor modifications, which are indicated as necessary. 

Theorem 4. Given a finite graph G = {V,E) andv e V, let N{v) = {vi, . . . ,Vk}- 
Then 

where G; is the graph induced by V \ {u, Ui, . . . , Vi}, Go is induced by V \ {v} 
and Jli<i<fe ~ 1 when fc = 0. 

Proof. We have 

I:v^I I:vel I:ICV\{v} I:ICV\{v,Vi,...,Vk} 

where everywhere the sums are over independent sets /. Note that Xl/-7cV'\-f«> ^ 
Zc„ and 
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J2i:icv\{v.vi,....v,,} ^'^' = ^Gfc- Dividing both sides of the identity (29l by Z^g 
we obtain 

It is immediate that Zq^/Zq = Pg(v ^ I)- In order to interpret Zq^JZiq^ 
similarly we rewrite it as 



-^(G 
^ 



n 



Zr, 



and note that Zq./Ziq._-^ 
obtain (pSl). 



j(wi ^ I). Combining these observations we 

n 



The identity ( 28 1 suggests a recursion for computing marginal probabilities 



{v ^ I) approximately. The idea is to apply the identity recursively several 
times and then set the initial values arbitrarily. One then establishes further 
that a correlation decay property holds on this recursion which implies that any 
initialization of the values at the beginning of the recursion leads to approxi- 
mately correct values at the end of the recursion. This principle was underlying 
the algorithm proposed in |Wei06] for computing approximately the number of 
independent sets in general graphs. The original approach taken in |Wei06] was 
slightly different - first a self-avoiding tree corresponding to the recursive com- 
putation tree described above is constructed. Then it is shown that P(w ^ /) 
on this tree equals the same probability in the underlying graph. The approach 
proposed here is slightly simpler as it bypasses the extra argument of showing 
equivalence of two marginal probabilities. 

We now provide details of this approach. Given a finite graph G, for every 
subgraph G = (V", E) of G, every vertex v €V and every i € Z_|_ we introduce 
a quantity <i>jj(w,t) defined inductively as follows. 



^g(^'*) 



1, 
(1 



A)- 



(i+ni<,<fe%_,(«»,i-i))~ 



t = 0; 




t > 0, N{v) = 0, 




t> 0, N{v) = {vi,. 


.,Vk}^^ 




(30) 



Here again Go is induced by T^\{w} and Gi is induced by V\{v, vi, . . . , Vi\. The 
recursion ( 30 1 is naturally related to the identity ( 28 1 . Specifically if <i>g {vi,t— 
1) = F^,_[vi ^ 7) for all i then ^^{v,t) = P^iyf. I)- However, this will not 
occur in general, as we set ^q{v, 0) = 1, due to the lack of knowledge of actual 
values of the corresponding probabilities. 

Similarly to <&, we introduce values \E'(j(?;,i) with the only exception that 
^(j(w, 0) = for all G — {V, E) and v ^V . The following lemma follows from 
Theorem |4] and the definitions of <i> and ^ using a simple induction argument 
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Lemma 1. For every v and t G Z_|_ 

*G(t',2t)<PG(«^I)<$G(«,2t), 

$g(«, 2t+l)< ¥g{v il)< *g(w, 2t + 1) 

Next we provide bounds on the computational effort required to compute <& 
and *. 

Lemma 2. For every finite graph G with degree < A, w G G and t, the values 
4>G(f , t), ^g(''^j ^) (^OLi^ be computed in time exp(0(ilog A)), where the constant 
in O(-) is universal. 

Proof. The result follows immediately from the recursive definitions of <& and 
^. D 

The crucial correlation decay property is formulated in the following propo- 
sition. 

Theorem 5 f [Wei06p . For every A > 3 and for every 

A<(A-l)^-V(A-2)^, (31) 

there exists C — C(A,A),p = /3(A,A) < 1 such that for every finite graph 
G = (y, -E) with degree at most A and every u G V, t G TLj^ : 

\log^G{v,t) ~log^G{v,t)\ < CpK (32) 



As a result, G satisfies exponential SSM for A satisfying (31) 



Proof. The details of the proof can be found in |Wei06| . It is shown there 



that the absolute value in ( 32 1 is upper bounded by the same quantity, when 
applied to G = Ta,* - the A-regular depth-i tree. Then a classical results by 
Spitzer |Spi75| and Kelly JKel85] are invoked to show the existence of C and 
p. The existence of p then implies SSM via Proposition [T] and observing that 
changing values of a for nodes u which have distance bigger than t from a given 
node w, does not affect the values of <I>G(f,t) and ^G(''^ji)- Q 

4.2 Free energy and surface pressure on Z'^. Numerical 
results 

We are now equipped to obtain bounds on the free energy and the surface 
pressure for the hard-core model on Z''. 

Denote by <I>(t) and ^(t) the values of <i>G(w,i), ^G(^',i) when applied to a 
graph G = Z^q n i?„, for sufficiently large n and v — Q. Observe that the values 
<I>g(i', t),^G(''^j^) E^re the same for all values of n sufficiently larger than t (for 
example n>t + l suffices) . Thus the notations are well-defined. The following 
relations are the basis for computing bounds on the free energy and surface 
pressure. 
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Corollary 3. For every t G Z+ and X satisfying iSl 



- log $(2t) < V{d, A) < - log *(2i) (33) 

'log^(2i+l) <rid,X) < -log$(2t+l). (34) 



Proof. By Theorem Is] we have exponential SSM. Thus 

(0^/)= lim P« nzi„(0^/)- 



^^0 



Bn n Z^Q is finite graph for which bounds from Lemma 111 are applicable. D 
Our algorithm for computing V{d, A) and s'P{'Z'^, A) is based on relations 



(33) and (34 1, and will be called Sequential Cavity Algorithm or shortly SCA. 

We have numerically computed values <i>(t),^(t) for the cases d = 2,3,4 
using the chessboard pattern method. Our results provide the following bounds 
on the free energy for the case A = 1. Since previous bounds were stated in 
terms of exp{V{d, 1), we do the same here: 

1.503034 < exp(P(2, 1)) < 1.503058 
1.434493 < exp(P(3, 1)) < 1.449698 
1.417583 < exp(-p(4, 1)) < 1.444713 

The computations were done at the level t = 27 for the case d — 2, t — 16 for the 
case d — 3 and t — 12 for the case d — i. Our lower bound for the case d = 2 is 
weaker than the previous best known 1.503047782 |CW98] . which is already very 
close to the presumably correct but unproven value stated in [Bax99 . However, 
our upper bound is stronger than the previous best known 1.5035148 |CW98] . 
We are not aware of any estimates for the case d = 3. Thus we believe our 
bounds are the best known. We have not done computations of the surface 
pressure and we are not aware of any previously existing benchmarks. Note 
that in the case d = 3 we have A = 6, and A = 1 no longer satisfies ( pl| . Thus 
we have no guarantee that SCA will provide converging estimates as t increases. 
The correctness of our bounds for this case is guaranteed by Corollary [l] It is 
encouraging to see that the bounds are close and based on this fact we conjecture 
that A = 1 corresponds to the uniqueness regime. 

It is instructive to compare our numerical results, which were obtained using 
the chess-pattern approach (identity (p7|) with results which could be obtained 



directly from (26 1. The computations based on ( |26[ ) for the case d = 2, A = 1 at 
depth t = 12 lead to bounds 1.0942 < exp('P(2, 1)) < 1.8377. At the same time, 
the computations using chess pattern method at the depth only t — 3 already 
lead to a much tighter bounds 1.4169 < exp('P(2, 1)) < 1.5565. 

Notice, that while Theorem [5] is not used in computing actual bounds on 
the free energy and surface pressure, it provides the guarantee for the quality 
of such bounds. Let us use it now to analyze the computation effort required to 
obtain a particular level of accuracy in bounds. 
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Proposition 2. For every d,X < (d— l)'^~^/((i— 2)'^ and e > SCA produces an 
e-additive estimate of'P{d, A) and sV{Z^, A) in time {-)^^\ where the constant 
in O(-) depends on A and d. 

Proof. Applying Theorem^ an additive error e is achieved provided that Cp* < 
e or 

t > log(C/e) log^^(l/p) — 0(log(l/e)). Hence the result for free energy follows 
from Lemma 121 

For surface pressure observe that applying Theorem [5] 

logP2;.^(ao = 0) - log P^. ^^^^((70 = 0) < CpK 

A similar bound holds for P^d (uo = 0). Thus if we take /cq such that 

Cp^° /{\ — p) < e, then the partial sum in ([7| corresponding to terms fc > fco is at 
most2e. The required fco is 0(log((l-p)C/e)/log(l/p)) = 0(log(l/e)). For the 
remaining terms fc < fcg we compute F^d ((Tq = 0) and P^d (ctq = 0) 

using SCA with accuracy e — e/fco- Since log(l/e) = log(l/e) + loglog(l/e) = 
0(log(l/e)), the result then follows from our estimate for computing Vid, A). D 

5 Monomer-dimer (matching) model 

The monomer-dimer model is defined by spin values S — {0, 1} assigned to 
edges of a graph G = {V,E). A set of edges M C -E is a matching if no two 
edges in M are incident. Sometimes term partial matching is used to contrast 
with full matching, which is a matching with size |T^|/2 (namely every node is 
incident to an edge in the matching). The edges of M are called dimers and 
nodes in V which are not incident to any edge in M are called monomers. We 
set iJ(0, 0) = H{0, 1) = H{1,0) = h{0) = 0, H{1, 1) = oo, h{l) = (i. The Gibbs 
measure is defined via (pi . As in the case of hard-core model, it is convenient 
to introduce A — exp(— /?) > 0. Similarly to the hard-core model we have for 
finite graphs G 



Zg = ^aI*^I (35) 



M 



were the summation is over all matchings in G. The summation in (35) is called 
a matching polynomial in the combinatorics literature. We denote by M a 
random matching chosen according to the Gibbs measure, when it is unique. In 
the case of finite graphs 

Pg(M = M) = Zg'aI*^I. 

The monomer-dimer model is a close relative of the hard-core model, even 
though its properties are substantially different. For example this model does 
not exhibit a phase transition and is always in the uniqueness regime [H L72] . 
Moreover, it satisfies the SSM property for all activities A as we shall shortly see. 
Corresponding analogues of Proposition [l] will be stated later once exponential 
SSM is asserted. 
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5.1 Recursion, sequential cavity algorithm and correlation 
decay 



We now establish an analogue of (|28|) for the monomer-dimer model. The proof 
of this result can be found in BGK+07 and is omitted. It is similar to the 



proof of (28 1. In the following, with a slight abuse of notation we write v G M 



if matching M contains an edge incident to v. 

Theorem 6. IBGK^ 01 ' For every finite graph G = {V, E) and v €V 



where Gq is induced by V \ {v} and X^neTVrftj) ~ ^ when v is an isolated node. 



The further development in this subsection mirrors the one of Subsection |4.1[ 
yet the conclusion will be different - the monomer-dimer model exhibits SSM 
for all values of A. This will lead to an algorithm for computing the free energy 
and surface pressure for monomer-dimer model for every value A > 0. 

Given a finite graph G, for every subgraph G = {V, E) of G, every node 
V e G and every t g Z_|_ we introduce a quantity <i>jj(w,t) defined inductively 
as follows. In the context of monomer-dimer model this quantity stands for 
(approximate) probability that v ^ M in the subgraph G. 

r 1, i = or Ar(t;) = 0; 

^gI^'^J- j (l + AEi<,<fc'i>Go(«.,t-l))^\ i>0and7V(^;) = {„^,...,„,}^ 

(37) 

Here Go is induced by 1/ \ {v}. If <i>Go(v»:* - 1) = Pc,„{vi ^ M) for all i then 
^G{v,t) = ¥G{v^M). 

Similarly, introduce ^(^(w, t) with the only exception that ^q(w, 0) = for all 
G and v. The following proposition follows from Theorem [6| and the definitions 
of <f> and ^ using a simple induction argument. 

Lemma 3. For every v E V, t e Z+ 

*G(«,2i) < Pg(v iM)< <i>G{v,2t) 
$g(w, 2t -f 1) < Pg(w iM)< ^g{v, 2t + 1). 

Next we provide bounds on the computational effort required to compute <I> 
and *. 

Lemma 4. For every finite graph G with degree A, w e G and t, the values 
4>G(f , t),^G('yj^) CO'''^ be computed in time exp(0(ilog A)), where the constant 
in 0{-) is universal. 

Proof. The result follows immediately from the recursive definitions of <f> and 
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The correlation decay property is formulated in the following proposition 
which is proved in |BGK+07| . Let 

/ 2 \ 1/2 

P = f 1 - , ) (38) 

Theorem 7 ([BGK^OT]). For every A > 2, A > 0, for every graph G with 
degree at most A and every node v 

|log$G(t',<)-log*G(w,i)|<P*log(l + AA). (39) 

As a consequence, every graph G satisfies exponential SSM for all A > 0. 

We now state the analogues of Proposition [T] and Corollary [2] The proofs 
are very similar and omitted. Given any set W C E, let N{W) — WU{e : 3e' G 
W,e-e'}. 

Proposition 3. Given a graph G = (V,E), for every mutually exclusive sets 
A,Wi,W2 C E and spin assignment {se),e G A on A the following identity 
holds with respect to the unique Gibbs measure. 

PG(l{e e M} == se, Ve £ A\Wi n M = 0, W^2 C M) = 

= Pe(l{eeM} = Se,Vee A), (40) 

where G is the subgraph obtained from G by removing edges Wi U N(W2) ■ 

We use notations V{d,X) and sV{d,\) for the free energy and the surface 
pressure for the monomer-dimer model on Z'' as well. As a corollary of Theo- 
rem [3] and Proposition [3] we obtain 

Corollary 4. For every A > and d 

V{d,\) = -logP^.^(0 iM) = -^logP^.^_ (0 i M) (41) 

5.2 Free energy and surface pressure on 7/. Numerical 
results 

We now obtain bounds on the free energy and surface pressure for the monomer- 
dimer model on lA. Again denote by <i>(i) and ^(i) the values of <i>G(f , t), ^g(i!, t) 
when applied to any graph G = Z^q n i?„ in the monomer-dimer context, for 



sufficiently large n and v ^Q. The relations ( |33| and (34 1 hold as well and the 
proof is very similar. 

Our algorithm for computing V(d, A) and sV{lA^ A) is again based on com- 
puting $(i) and \E'(i) and is again called Sequential Cavity Algorithm (SCA). We 
now report numerical results on computing 'P{d, A). We have computed values 
4>(i),^(i) for d— 2,3,4 and a range of values A using the chessboard pattern 
method. Our upper and lower bounds for free energy are presented in Table [l] 
The depth levels t — 14,9,7 were used for the cases d — 2,3,4. As expected. 
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our bounds are high quahty for lower A and then degrade as A — + oo. Each 
computation run took about 3 minutes on a workstation and we have not made 
an attempt to obtain very accurate bounds for each value A. However for the 
case of interest A = 1 we ran our algorithm for larger depths. For the case d=2 
a very accurate rigorous estimate 0.6627989727 ± 0.0000000001 is due to Fried- 
land and Peled |FP05] . shown non-rigorously earlier by Baxter |Bax68] . We 
have not made an attempt to improve this bound. However for the case d = 3 
we can significantly improve the best known bound 0.7850 < 7^(3, 1) < 0.7863 
due to Friedland et al |FP05j . |FG08] . |FKLM] . At depth t = 19 we obtained 
estimates 0.78595 < 7^(3, 1) < 0.78599, which is two orders of magnitude im- 
provement. One should also note that the lower bounds 0.7845 in |FG08] and 
0.7850 in |FKLM| were obtained using Friedland-Tveberg inequality which pro- 
vides a bound for general regular graphs. Thus, while highly accurate for the 
case Z^, this bound is not improvable by running some numerical procedure 
longer or on a faster machine. The previous best known numerical estimate 
7^(3, 1) > 0.7653 [FP05^, is based on the transfer matrix method is weaker. We 
have also obtained bounds for d = A for which no prior computations are avail- 
able. We obtained 0.8797 < 7^(4, 1) < 0.8812. The computations were done at 
depth t = 14. 

The following proposition gives a bound on the numerical complexity of 
SCA. 

Proposition 4. For every d > 2,X > 0,e > SCA produces an e-additive 

estimate of'P{d, X) and s'P{d, X,a) in time log{l + 2Xd){l/ e) ^ ' °® ' , where 
the constant in 0{-) is universal. 

Again for constant X,d we obtain performance (l/e)'^(^\ which is a quali- 
tative improvement over the numerical effort exp(0((l/e)''~^)) of the transfer 
matrix method. 

Proof. We have A = 2d. Applying Theorem |7] an additive error e is achieved 
provided that p* log(l + 2Xd) < e or 



t> (-logp)"i(log-+loglog(l + 2Ad)) 



Applying g we have logp = 0( ^i+L+i ) = 0{^;). The result for P(d, A) 
then foUowsirom this estimate and Lemma H The result for sV{d, A) is estab- 
lished using the same line of reasoning as for Proposition |2] D 

6 Conclusions 

Several statistical mechanics models besides hard-core and monomcr-dimer mod- 
els fit our framework, yet were not discussed in this paper. One such model is 
Ising model and its generalization. Potts (coloring) model |Sim93] . The Ising 
model satisfies our Assumption [T] by making s* = 1 or —1. Thus the represen- 
tation Theorem fl] holds, provided (exponential) SSM holds. However, in order 
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A 


d=2L 


d = 2V 


d = 3L 


d=3U 


rf = 4L 


(i:=4U 


0.1000 


2.3219 


2.3219 


2.3311 


2.3311 


2.3399 


2.3399 


0.2000 


1.6802 


1.6802 


1.7096 


1.7096 


1.7363 


1.7363 


0.3000 


1.3451 


1.3451 


1.3959 


1.3959 


1.4395 


1.4395 


0.4000 


1.1361 


1.1361 


1.2050 


1.2050 


1.2621 


1.2624 


0.5000 


0.9934 


0.9934 


1.0770 


1.0770 


1.1442 


1.1453 


0.6000 


0.8902 


0.8902 


0.9853 


0.9855 


1.0599 


1.0629 


0.7000 


0.8122 


0.8122 


0.9164 


0.9171 


0.9961 


1.0025 


0.8000 


0.7513 


0.7513 


0.8627 


0.8642 


0.9457 


0.9568 


0.9000 


0.7026 


0.7026 


0.8196 


0.8224 


0.9044 


0.9215 


1.0000 


0.6628 


0.6628 


0.7840 


0.7887 


0.8695 


0.8937 


1.5000 


0.5389 


0.5392 


0.6677 


0.6890 


0.7475 


0.8163 


2.0000 


0.4742 


0.4758 


0.5982 


0.6436 


0.6671 


0.7840 


3.0000 


0.4056 


0.4139 


0.5079 


0.6055 


0.5564 


0.7586 


4.0000 


0.3665 


0.3855 


0.4455 


0.5908 


0.4793 


0.7492 


5.0000 


0.3390 


0.3704 


0.3972 


0.5836 


0.4208 


0.7448 


6.0000 


0.3174 


0.3614 


0.3575 


0.5797 


0.3740 


0.7423 


7.0000 


0.2993 


0.3558 


0.3239 


0.5772 


0.3354 


0.7408 


8.0000 


0.2837 


0.3519 


0.2947 


0.5757 


0.3026 


0.7399 


9.0000 


0.2697 


0.3493 


0.2691 


0.5746 


0.2743 


0.7392 


10.0000 


0.2570 


0.3473 


0.2465 


0.5738 


0.2496 


0.7387 


15.0000 


0.2064 


0.3426 


0.1643 


0.5719 


0.1621 


0.7376 


20.0000 


0.1685 


0.3409 


0.1147 


0.5713 


0.1110 


0.7372 


30.0000 


0.1151 


0.3397 


0.0627 


0.5708 


0.0592 


0.7369 


40.0000 


0.0811 


0.3393 


0.0386 


0.5706 


0.0359 


0.7368 


50.0000 


0.0592 


0.3391 


0.0259 


0.5705 


0.0239 


0.7368 



Table 1: Upper (U) and lower (L) bounds on free energy for d = 2,3,4 
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to turn it into a useful method for obtaining provably converging bounds, we 
need an analogue of Theorems Is] and IT] namely the correlation decay property 
on a computation tree. Such a result is indeed established in [GK07] . but for 
an Ising model with very weak interactions. Thus it is of interest to strengthen 
the result in |GK07] and obtain some concrete estimates for the Ising model 
on 1^. We are not aware of earlier benchmark results for this model. The 
situation with the Potts model is similar with the exception of hard-core Potts 
model (proper coloring). In this case the required s* does not exist, as for every 
color s, -ff (s, s) = oo. However, our method can be extended by considering for 
example a periodic coloring of nodes v >- ^ ox using the chess-pattern version 
Theorem [2] of our main result. The required correlation decay result is estab- 
lished in |GK07] for the case q > 2.86A and q, A appropriately large constants. 
Thus in order to turn this result into a useful method for computing free energy 
and surface pressure for coloring model on Z'^, one needs to deal with these con- 
stants explicitly. Also a convenient monotonicity present in the hard-core and 
monomer-dimer models is lost in the Potts model case, which makes application 
of Corollary [1] harder. 
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